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ABSTRACT 


This  thesis  addresses  three  deficiencies  in  the  NPSNET  simulated  world.  First, 
although  a  full  set  of  algorithms  have  been  defined  for  Dead  Reckoning  (DR)  entities  in  a 
distributed  simulation,  NPSNET  only  implements  a  few  simple  linear  algorithms.  Second, 
NPSNET  lacks  a  set  of  physically-based  models  for  munition  trajectories  having,  currently, 
only  algorithms  for  the  bullet  and  bomb.  Third,  NPSNET  lacks  physically-based  models 
for  engine  power  curves  using,  instead,  a  simple  linear  approximation. 

The  purpose  of  this  thesis  work  is  to  implement  an  object-oriented  programming 
toolkit  which  corrects  these  deficiencies.  The  code,  in  C++,  utilizes  class  hierarchies.  The 
toolkit  implements  the  nine  class  hierarchies  of  DR  algorithms  described  by  the  Advanced 
Research  Projects  Agency  for  the  Distributed  Interactive  Simulation  standard.  The  toolkit 
also  provides  treatment  of  a  physically-based  class  hierarchy  for  munitions  trajectories.  In 
addition,  a  physically-based,  mathematical  model  for  the  engines  class  was  implemented. 

As  a  result,  a  set  of  DR  algorithms  have  been  built  to  predict  the  position  of 
simulated  entities  in  all  cases.  The  munitions  class  implements  trajectories  for  a  variety  of 
projectiles.  With  this  arsenal,  future  versions  of  NPSNET  will  be  more  realistic.  The  engine 
class,  with  new  mathematical  models,  far  more  realistically  represents  engine  behaviors 
than  the  current  linear  approximation.  In  sununation,  the  implementation  of  this  toolkit 
dovetails  very  well  with  the  needs  of  NPSNET,  and  wiU  be  integrated  into  future  releases. 
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1.  INTRODUCTION 


A.  BACKGROUND 

Conventional  training  and  preparation  for  future  military  operations  is  an  expensive 
and  time-consuming  process  involving  many  persons  and  resoiuces.  With  the  current 
budgetary  constraints  and  downsizing  of  the  armed  forces,  ways  of  doing  business  must  be 
revised,  but  in  a  way  that  does  not  compromise  preparedness.  Fortunately,  with  advances 
in  Computer  Science,  military  operations  from  battle  to  replenishment  can  now  be 
simulated  on  computers,  and  personnel  can  complete  many  phases  of  training  in  simulated 
battlefields,  thus  mitigating  the  expense  of  engaging  in  the  real-world  environment. 

One  such  simulated  battlefield,  NPSNET  [PRAT93],  has  been  developed  by  the 
Computer  Science  Department  of  the  Naval  Postgraduate  School  (NPS).  NPSNET  is  a  low- 
cost,  distributed,  interactive,  virtual  3D  world  for  simulation  and  training.  Personnel  in  the 
simulated  environment  interact  with  each  other  as  if  they  were  on  a  battlefield.  They  can 
virtually  fly  an  airplane,  ride  a  tank,  walk,  run,  and  shoot  at  each  other. 

NPSNET  is  an  evolving  system,  a  test-bed  for  new  research  as  well  as  a  growing 
amalgam  of  proven  ideas.  Recently,  numerous  projects  have  improved  and  enhanced 
NPSNET.  Such  works  encompass:  Configuration  Management,  Real-time  Scene 
Management,  Terrain  Development,  Insertion  of  the  Human,  Human-Computer  Interface, 
and  Networking.  However,  NPSNET  lacks  trae  physical  representation  of  moving  bodies. 
Only  a  few  simple  cases  have  been  implemented  such  as  linear  Dead  Reckoning  (DR),  and 
simplified  munitions  trajectories  for  the  bullet  and  bomb.  Improving  the  physically-based 
aspects  of  NPSNET  is  the  purpose  of  this  thesis  work. 

B.  PROBLEM  STATEMENT 

This  thesis  addresses  and  resolves  three  specific  deficiencies  in  the  physically- 
based  representations  of  NPSNET.  The  first  is  dead  reckoning  (DR).  NPSNET  uses  the 
Distributed  Interactive  Simulation  (DIS)  protocol,  a  standard  for  networked  simulators,  to 
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communicate  between  entities-aircraft,  vehicles,  munitions,  etc.  -  in  its  environment 
[IST93].  Every  workstation  participating  in  a  simulation,  each  running  NPSNET,  shares 
the  data  for  the  entities  it  controls,  thus  allowing  each  NPSNET  to  display  a  consistent  view 
of  the  simulated  world  to  its  user.  Each  entity  in  the  environment  sends  to  the  network  its 
velocity,  position,  and  orientation  no  less  frequently  than  every  five  seconds.  Other  objects 
receive  and  render  that  information  within  their  own  workstation.  For  efficiency  and  to 
conserve  network  bandwidth,  data  for  an  entity  is  not  sent  out  continually.  To  keep  an 
updated  picture  of  networked  entities,  each  workstation  dead  reckons  (predicts)  remote 
entities’  velocity,  position,  and  orientation  based  on  its  last  update.  Specific  dead  reckoning 
algorithms  constitute  part  of  the  DIS  standard.  Currently,  NPSNET  implements  only  a  few 
simple  linear  world  coordinate  cases.  A  full  implementation  of  DR  for  all  vehicle  cases  in 
both  world  and  body  coordinates  is  required. 

The  second  deficiency  regards  simulation  of  munitions.  In  a  simulated  battle, 
continuous  simulation  of  munitions  is  absolutely  essential.  Currently  in  NPSNET,  only  a 
few  simple  cases  of  a  bullet  and  a  bomb  are  implemented.  A  more  definitive  physically- 
based  implementation  of  munitions  is  required  to  more  closely  reflect  real-world  munitions 
behavior. 

The  third  deficiency  is  engine  simulation.  To  make  the  virtual  world  realistic  and 
plausible  to  the  user,  moving  objects  must  reflect  the  behavior  of  their  real-world 
counterparts  as  much  as  possible.  Currently,  NPSNET  lacks  any  engine  type.  Therefore,  a 
basic  implementation  of  engines  is  needed,  too. 

C.  APPROACH 

For  DR  implementation,  we  use  the  study  of  the  Advanced  Research  Projects 
Agency  (ARPA)  for  DR  algorithms  published  as  part  of  DIS  2.0.3  [TOWE94].  In  this 
study,  the  equations  of  motion  for  the  DR  algorithms  are  divided  into  nine  groups.  The  first 
group  is  for  a  static  case.  Another  four  groups  are  for  world  coordinates,  and  the  remainder 
are  for  body  coordinates.  Using  the  C++  object-oriented  programming  methodology,  this 
thesis  work  implements  all  nine  groups  into  a  class  hierarchy. 
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The  munition  implementation  is  based  on  Newtonian  Mechanics.  Most  moving 
objects  in  the  real  world  are  governed  by  gravitational  forces.  Therefore,  kinematic 
equations  with  constant  downward  acceleration  are  utilized  for  bullet,  bomb,  ballistic,  and 
artillery  objects.  A  physically-based  model  of  a  system  with  lost  mass  is  used  for  rockets 
and  missiles.  Air  drag  is  omitted. 

For  engine  implementation,  we  built  an  abstract  input  and  output  system  with 
general  quadratic,  logarithmic,  and  exponential  models  of  acceleration.  (A  full  physically- 
based  model  is  beyond  the  scope  of  this  thesis.) 

All  of  the  implemented  code  is  written  in  C-h-  with  class  hierarchies  for  portability 
and  reusability.  With  its  class  structure  and  inheritance,  C++  minimizes  code  development 
as  well  as  making  the  program  easier  to  understand. 

D.  ORGANIZATION 

This  first  chapter  introduced  NPSNET  and  discussed  the  problems  addressed  by 
this  thesis.  It  also  described  the  approach  taken  to  solve  each  problem.The  second  chapter 
covers  the  background  and  previous  works  impinging  on  the  problems  to  be  solved.  The 
third  chapter  describes  in  detail  the  implementation  of  the  DR  standard.  The  fourth  chapter 
recounts  the  formulas  and  equations  which  govern  the  munition  classes.  The  fifth  chapter 
describes  the  implementation  of  the  engine  mathematical  models  and  the  derived  formulas. 
The  sixth  chapter  contains  the  conclusion  and  presents  necessary  future  work. 
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11.  BACKGROUND  AND  PREVIOUS  WORK 


A.  PREVIOUS  WORKS 

1.  Dead  Reckoning  (DR) 

Currently,  NPSNET  implements  only  a  few  simple  cases  of  DR  where  there  is  no 
rotation  involved.  There  are  the  static  case,  the  two  cases  of  world-coordinates,  and  the 
default  case  for  the  rest  of  the  situations,  see  Figure  1. 


Figure  1:  Current  NPSNET  DR  IMPLEMENTATION 


The  first  algorithm  applies  to  a  static  object  when  its  rotation,  velocity,  and 
acceleration  are  all  zero.  The  second  algorithm  applies  when  object  rotation  and 
acceleration  are  zero  but  velocity  is  not.  The  third  algorithm  applies  when  object  rotation 
is  zero  and  velocity  and  acceleration  are  not.  The  default  case  applies  when  object  velocity 
is  a  constant,  so  that  the  object  moves  in  a  straight  line.  These  four  DR  cases  were 
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implemented  in  C++  as  special  cases  in  the  subclass  moveDR  of  the  Vehicle  Class  in 
NPSNET.  Kinematic  equations  with  constant  acceleration  were  used  to  implement  those 
cases: 

V;  =  V.  +a;t  Eql 

I  t. 

X;  =  X;  +v,i  +  0.5nr^  Eq2 

^  '■o  ‘ 

Where: 

V-  =  respective  velocity  in  each  direction 
v-  =  respective  initial  velocity  in  each  direction 
jCj.  =  respective  position  in  each  direction 
X  =  respective  initial  position  in  each  direction 
a-  =  respective  acceleration  of  each  direction 
t  =  time 

One  of  the  purposes  of  DIS,  a  set  of  standards  developed  by  ARPA  and  industry,  is 
to  provide  a  specification  to  be  used  by  government  agencies  and  engineers  to  build 
interoperable  networked  simulation  systems  [IST93].  DIS  defined  a  need  for  DR,  and 
ARPA  released  its  DR  algorithms  study  on  February  7, 1994  [TOWE94].  In  this  study,  with 
C++  methodology  in  mind,  the  algorithms  are  divided  into  a  general  static  case  and  two 
large  groups,  one  for  the  world-coordinates,  and  the  other  for  body-coordinates.  Inside  each 
large  group  are  four  subgroups,  each  addressing  different  cases  of  the  moving  object  in  its 
respective  coordinates,  see  Figure  2. 
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World  coord. 


Body  coord. 


Figure  2:  DR  algorithms  of  ARPA’s  study 


2.  Munition 

In  the  current  version  of  NPSNET,  only  a  few  simple  cases  of  munitions  are 
implemented,  such  as  a  bullet,  a  bomb,  and  a  simple  non-loss-mass  missile,  see  Figure  3. 
Equations  similar  to  Eql  and  Eq2  above  are  used  to  implement  those  munitions. 


Figure  3:  Current  NPSNET  munitions  implementation 


Newtonian  Mechanics  governs  the  behavior  of  objects  flying  in  the  air  in  a  constant 
gravitational  field.  Equations  that  describe  motion  and  kinematics  have  been  developed 
ever  since  Newton  presented  his  greatest  discovery  (F  =  ma).  Many  first-year  math  and 
physics  books  describe  these  equations.  [MART84][HALL78] 


Another  relevant  work  is  the  NATO  study  of  ballistic  missiles  [NATO].  It  is  the 
definitive  study  of  the  trajectory  of  a  flying  mass  in  the  air.  It  describes  all  factors 
influencing  the  behavior  of  flying  objects,  encompassing  drag  force,  lift  force,  Magnus 
force,  and  Coriolis  force.  Most  of  these  equations  exceed  the  scope  of  this  thesis  and  are 
unnecessary  for  NPSNET,  though  it  could  apply  to  later  enhancement  work. 

3.  Engines 

Engine  models  are  not  currently  incorporated  into  NPSNET.  Presently,  NPSNET 
uses  a  simple  linear  approximation  between  the  user  input  to  a  throttle  and  the  speed  of  the 
entity.  However,  the  real  world  of  vehicles  encompasses  a  variety  of  engine  models,  such 
as  motors,  gas  turbine,  piston,  and  rocket  engines.  All  currently  existing  engine  models 
exceed  the  scope  of  this  thesis,  except  for  the  rocket  engine.  With  a  few  revisions,  it  follows 
the  rocket  munition  model  and  so  is  implemented  in  this  thesis  work. 

B.  SUMMARY 

In  summation,  NPSNET  lacks  complete  physically-based  models  for  DR, 
munitions,  and  engines,  for  accurately  representing  objects  in  its  virtual  world.  All  DR  and 
munition  models  are  readily  available,  and  can  be  implemented.  Except  for  the  rocket 
engine  which  can  be  implemented  directly  from  rocket  models,  the  remaining  engine 
models  need  to  be  abstracted  as  an  input  and  output  system  and  tested  until  it  seems 
plausible  in  the  virtual  world. 
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m.  DEAD  RECKONING 


A.  OVERVIEW 

What  is  dead  reckoning?  According  to  the  Random  House  dictionary,  dead 
reckoning  is  the  “calculation  of  one’s  position  on  the  basis  of  distance  run  on  various 
headings  since  the  last  precisely  observed  position.”  Our  usage  of  the  word  dead  reckoning 
in  this  chapter  is  exactly  the  same  as  the  definition. 

To  allow  multiple  users  in  a  simulation,  we  ran  NFS  NET  on  different  workstations 
which  are  aU  connected  by  the  same  network.  Each  workstation  presents  a  consistent  view 
of  the  entire  simulated  world,  rendering  both  the  local  objects  it  controls  and  the  networked 
objects  controlled  by  the  other  simulators.  Each  simulator  transmits  the  position,  velocity, 
and  orientation  for  its  own  local  objects,  and  reads  the  network  for  information  on  remote 
objects.  This  state  change  data  is  sent  on  the  network  no  less  frequently  than  every  five 
seconds  in  a  packet  called,  in  DIS  terminology,  a  Program  Data  Unit  (PDU).  But  five 
seconds  is  a  long  time  when  the  user  is  watching  an  entity  which  is  animated  on  the  screen. 
The  entity  cannot  jump  to  a  new  location  every  so  many  seconds  and  satisfy  the  appearance 
of  smooth  motion.  In-between  times,  each  workstation  needs  to  dead  reckon  each  remote 
object’s  position,  velocity,  and  orientation  and  move  it  smoothly  along  its  probable  course. 
Figure  4. 


PDU 


DR  DR  DR 

- 5  seconds - 


DR 


PDU 


Figure  4:  PDU  and  DR  time 
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In  endeavoring  to  provide  a  complete  set  of  DR  functionality,  this  thesis  utilizes  the 
algorithms  developed  by  ARPA  for  both  world  and  body  coordinates.  These  are 
implemented  using  appropriate  C.++  classes. 


Figure  5:  DR  classes 


In  object-oriented  languages  like  C++,  inheritance  is  when  a  child  class  inherits  the 
variables  and  functions  of  a  parent  class.  Thus  the  child  class  is  automatically  like  its  parent 
except  for  the  specific  ways  in  which  it  needs  to  be  different,  which  makes  design  and 
coding  clean  and  efficient.  In  Figure  5,  each  of  the  nine  ellipses  represents  a  class  of  DR 
algorithm.  The  parent  class  is  the  static  class.  The  four  classes  on  the  left  side  of  the 
structure  belong  to  the  world  coordinates  class,  the  four  on  the  right  belong  to  the  body 
coordinates. 

The  structure  in  Figure  5  varies  a  bit  from  ARPA’s  DR  study,  see  Figure  2.  As 
represented  in  Figure  2,  dead  reckoning  method  four  (DRM4)  can  be  derived  firom  either 
DRM3  or  DRM5;  likewise  for  DRM8,  which  can  be  derived  from  DRM7  or  DRM9.  But, 
in  Figure  5,  DRM4  is  only  derived  firom  DRM3,  and  DRM8  is  only  derived  from  DRM7. 
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The  reason  for  the  variance  is  to  reduce  the  complications  arising  from  multiple  inherited 
variables  and  functions.  It  is  far  more  straightforward  to  keep  track  of  all  the  variables  and 
functions  if  each  child  class  derives  from  a  single  parent  class. 

Table  1  lists  each  DR  algorithm  in  the  set.  A  key  to  the  abbreviations: 

B  =  body  coordinate 
F  =  fixed 
P  =  position 
R  =  rotation 
V  =  velocity 
W  =  word  coordinate 

For  example,  the  fourth  row  indicates  that  algorithm  DRM4  has  rotation,  velocity, 
and  belongs  to  the  world  coordinates  subclass.  Its  roll,  pitch,  yaw  (orientation),  velocity, 
and  acceleration  are  equal  to  some  number. 


n 

H 

i 

■fS' 

1 

Static 

0 

0 

0 

0 

0 

2 

DRM2(F,P,W) 

0 

0 

0 

V 

0 

n 

:> 

DRM3(R,P,W) 

roll 

pitch 

yaw 

V 

0 

4 

pitch 

V 

a 

5 

DRM5(F,V,W) 

0 

0 

V 

a 

6 

DRM6(Fd>,B) 

0 

0 

0 

V 

0 

7 

DRM7(R,P,B) 

roll 

pitch 

yaw 

V 

0 

8 

DRM8(R,V,B) 

roll 

pitch 

yaw 

V 

a 

9 

DRM9(R,V,B) 

0 

0 

0 

V 

a 

Table  1.  DR  Algorithms 


Furthermore,  in  the  representation  DRM(X,Y,Z),  X  represents  for  the  rotation  of 
the  entity,  Y  represent  for  position,  and  Z  represent  for  coordinate. 
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C.  IMPLEMENTATION 


I.  Static  class  (DRMI) 

This  is  the  parent  class  for  all  other  DR  classes.  Its  rotation,  velocity,  and 
acceleration  are  all  zero,  so  its  position  is  fixed.  To  use  it,  we  simply  instantiate  an  object 
and  set  its  position  to  given  values.  In  this  case,  we  always  get  back  the  same  initial  value, 
when  we  call  for  updated  position. 


2.  World-coordinates  (DRM2-5)  classes 
The  input  to  all  the  algorithms  are: 

j  =  position  vector  in  world  coordinates  at  initial  time, 
j  =  velocity  vector  in  world  coordinates  at  initial  time. 

=  acceleration  vector  in  world  coordinates. 

( V  ( ,0  ( ( ^o) )  =  Euler  angles  at  initial  time. 

(0  =  ^(0^,0)^, =  angular  velocity  in  body  coordinates. 
A/  =  time  increment  for  dead  reckoning  step. 


The  outputs  are: 

^w[^o  j  “  position  vector  in  world  coordinates  after  time  delta. 
+  Atj  =  velocity  vector  in  world  coordinates  after  time  delta. 
+  Arj,(j)(^r^  +  Atjj  =  Euler  angles  after  tir 


time  delta. 


(1)  DRM2 

In  this  algorithm,  r  =  0,  v  =  x,  and  a  =  0.  We  use  the  Eq3  and  Eq4  to  compute 
the  object  velocity  and  position  at  time  t^  +  At.  Euler  angles  at  time  delta 

+ Atj,0^t^  + Atj,(t)^r^  + Arjj  are  the  same  as  its  initial  time 

jj  because  there  is  no  rotation  involved. 


12 


V  It  +An  =  7  \t  1  +  AM 
w\  o  J  w\  o  ]  w 


X  it  +At]=X  t  ]+AtV  it  +^A?X, 
w\  o  )  w\  o  ]  w[  o  2  w 


(2)  DRM3 


In  this  algorithm,  r  =  x,  v  =  x,  and  a  =  0.  We  can  get  velocity  and  position  at 
time  t^  +  At  by  Eq3  and  Eq4,  and  the  steps  to  get  Euler  angles  at  time  t^  +  At  are  as 

follows: 

Step  1:  Use  the  initial  Euler  angles  (\|/  (t^)  ,6  (i^)  (t^) )  to  compute  the  rotations 

matrix  U  (b  |  w)  which  takes  world  coordinates  into  body  coordinates  at  time  t^ .  The 
formula  for  U  (b  (tj  |w)  is 

c(\|f)c(0)  i(\|f)c(e)  -s(e) 

H  A ~  +s(<t')«(6)c(v)  c((t))c(v)  +s((l))s(e)s(v)  s(<t))c(e)  EqS 

[5((j))s(\|f) +c((!))s(e)c(\i/)  — s((l))c(v) +c((l»)c(0)s(\|/)  c((i))c(0) 

Where  c  =  cosine  and  s  =  sine. 

Step  2:  Compute  the  rotation  matrix  U  (b  (t^  +  At)  |b  (tj )  which  takes  body 
coordinates  at  time  t^  into  body  coordinates  at  time  t^  +  At. 

u(b(,  tA.'lli.f,  11  =  ( 1  -  cos  (N  AO)  ^  J  ^  .  sil.  (Iml  A.)  „  Eq6 


Where: 


0  -CO3  0)2 


=  CO3  0  -co^ 


-0)2  o)j  0 


The  matrix  I  is  the  3x3  identity  matrix  having  I’s  on  the  main  diagonal  and  O’s  off 


of  the  main  diagonal.  The  matrix  coto^  is 
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CO^COj  CDjCO^  COjCOg 


coco^  = 


“2®2  “2®3 
0)30)^  ©3^2  *^3^3 


and  the  magnitude  of  co  is 

|co|  =  (co^  +  (o|  +  co|) 


Eq8 


Eq9 


Step  3:  Compute  the  rotation  matrix  U  {b(t^  +  At)^w)  which  takes  world 
coordinates  into  body  coordinates  at  time  t^  +  At. 

Uibit^  +  At)\w)  =  U(bit^  +  At)\bit^))  ■UibUo)\^)  EqlO 

Step  4:  Compute  the  Euler  angles  (\|/i^  +  Ai)  +  At)  ,(|)  (t^  +  At) )  after  the 
time  increment  using  the  rotation  matrix  U  (b  (t^  + At)  ^w) .  The  rotation  matrix 
U  (b  (t^  +  At)  ^w)  has  the  form; 

U{b(t^  +  At)\w)  = 

Where  the  entries  U^j  are  functions  of  the  Euler  angles 
(\|/  +  Ai)  ,6  shown  in  the  matrix  displayed  in  Step  3  above. 

Step  5:  The  Euler  angle  0  is  recovered  via 

e  =  -asin(£/i3)  e  Eql2 

If  and  do  not  both  equal  zero,  the  Euler  angle  \|/  is  recovered  via 
\|r  =  Arg{U^.^,U^j)  e  [-Jt,7c]  Eql3 


^11  ^12  ^13 
^21  ^22  ^23 
C/31  ^32  ^33 


Eqll 
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Where  the  function  Arg(x,y)  is  an  argument  of  the  complex  number  x  +  yi,  the 
particular  branch  of  the  argument  chosen  to  lie  in  the  interval  [-7C,7C]  . 

If  and  £7^^  do  not  both  vanish,  the  Euler  angle  ^  is  recovered  via 

(j)  =  Arg(  1/33, t/23)  Eql4 

The  exceptional  cases  occur  when  6  =  7t/2  or  -ic/2 .  In  that  case,  the  Euler  angles 
\|r  and  ^  are  recovered  via 

\|/  =  0  Eql5 

(])  =  Arg(£/22“^32) 

(1)  DRM4 

In  this  algorithm,  r  =  x,  v  =  x,  and  a  =  x.  Position,  velocity,  and  orientation 
are  computed  much  the  same  as  DRM3  but  with  the  addition  of  acceleration  A^ . 

(3)  DRM5 

In  this  algorithm,  r  =  0,  v  =  x,  and  a  =  x.  We  can  use  Eq3  and  Eq4  above  to 
compute  its  position  and  velocity.  DRM5  lacks  rotation  so  we  do  not  require  equations  Eq5 
through  Eql6.  Its  orientation  is  fixed. 

3.  Body-coordinates  (DRM6-9)  classes 

The  input  to  all  the  algorithms  are: 

(t^)  =  position  vector  in  world  coordinates  at  initial  time. 

=  velocity  vector  in  body  coordinates  at  initial  time. 

A^  =  time  derivative  of 

( V  ( ^o)  ,0  ( ^(,)  =  Euler  angles  at  initial  time. 

(0  =  ^ CO  j, 002,0)3  =  angular  velocity  in  body  coordinates. 

At  =  time  increment  for  dead  reckoning  step. 
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The  outputs  are: 


^wyo  ^  ~  position  vector  in  world  coordinates  after  time  delta. 

+  At  j  =  velocity  vector  in  body  coordinates  after  time  delta. 

+  At  +  Atj,(t)^t^  +  At  j  j  =  Euler  angles  after  time  delta. 

(4)  DRM6 

In  this  algorithm,  r  =  0,  v  =  x,  and  a  =  0.  We  can  compute  the  velocity  and 
position  at  t^  +  At  with  Eql7  and  18,  ignoring  second  terms  on  the  right  of  both  equations 

because  they  are  zero.  The  Euler  angles  at  time  t^  +  At  is  the  same  as  its  initial  condition 

(¥  Uq)  >6  ,(})  (t^) ) ,  because  there  is  no  rotation  involved. 


ri<,xQ*  ,  |m|At  -  sin.(MMo)ar^  (IM ^0  ^  ^  1  -  cos  (NAQ  ^ 

JO  lfn|3  0)  0)2 


i|o)pAt^-cos  (|o)|At)  -lo)|Atsin  (|o)lAt)  +  1 


0)0)^  + 


cos  (|a)|At)  +  |o)|Atsin  (|a)|At)  - 1 ,  .  sin  (|(o|At)  -  |to|Atcos  (|o)|At)  ^  ^ 

M2  W 

Computing  Ml  and  il  is  the  same  as  in  Eq9  and  Eq7,  and  ttf  is  the 

transposition  of  the  matrix  U  (b  (t^)  jw)  above.  See  Eq5. 
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The  following  steps  are  used  to  calculate  Euler  angles 

(\|/  +  AO  ,e  +  At) ,(()  +  At) )  at  time  t^  +  At. 

Step  1;  Compute  the  rotation  matrix  U  (b  (t^  +  At)^(t^))  which  takes  body 
coordinates  at  time  t^  into  body  coordinates  at  time  t^  + At,  the  same  as  Eq  6. 

Step  2:  Compute  the  rotation  matrix  U  {b{t^  +  At)^^)  which  takes  world 
coordinates  into  body  coordinates  at  time  +  At  as  Eq  10. 

Step  3:  Compute  the  Euler  angles  (\|t  {t^  +  AO  ,6  +  AO  +  At) )  after  the 

time  increment  using  the  rotation  matrix  U  {b{t^  +  At)  jw)  ,  same  as  the  Step  4  in  the 

world  coordinates. 

(5)  DRM7 

In  this  algorithm,  r  =  x,  v  =  x,  a  =  0.  We  can  compute  the  velocity,  position, 
and  Euler  angles  at  t^  +  At  with  Eq3,  Eq4,  and  Steps  1  to  3. 

(1)  DRM8 

The  steps  to  compute  the  outputs  for  this  algorithm,  where  r  =  x,  v  =  x,  and 
a  =x,  are  the  same  as  the  steps  of  DRM7  with  the  exception  that  a  is  not  zero. 

(6)  DRM9 

In  this  algorithm,  r  =  0,  v  =x,  and  a  =x.  The  steps  to  get  the  outputs  are  the 
same  as  DRM6  with  a  not  equal  to  that  zero.  The  only  thing  changes  in  this  algorithm  is 
that  the  second  term  on  the  right  of  Eql7  is  not  zero. 

D.  RESULTS 

With  the  organization  afforded  by  the  C++  class  hierarchy,  the  coding  and  testing 
procedures  for  DR  classes  are  segregated  into  separate,  discrete  units.  The  debugging 
process  for  each  implemented  algorithm  was  very  simple.  For  each  case,  the  code  was 
tested,  as  suggested  in  Figure  4,  as  follows  (see  code  in  Appendix  A  for  details); 
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a.  Instantiate  an  object  of  that  case. 

b.  Call  procedure  update  DR  to  set  the  object  to  its  “real”  values  (as  the  information 
is  obtained  from  incoming  PDUs). 

c.  Intermittently  call  move  DR  to  compute  the  predicted  position,  velocity,  and 
angular  position  of  the  object. 

d.  Repeat  step  c.  at  the  desired  interval  until  the  next  PDU  arrives  with  real  data  for 
the  object,  then  repeat  step  b. 

The  code  is  extremely  reliable,  and  its  performance  very  quick,  even  when  tested 
on  a  low-end  Intel  80486-based  PC.  To  illustrate  test  results,  the  following  Tables  2  and  3 
list  the  input  and  output  tested  cases  for  algorithms  DRM4  and  DRM8.  They  are  the  two 
most  complex  cases  in  the  DR  classes. 


Initial  Conditions 

X 

y 

z 

Position  -  m 

-6378137 

-5 

-10 

Velocity-  m/sec 

30 

35 

40 

Acceleration  -  m/sec*sec 

-80 

-85 

-90 

Orientation  -  deg 

15 

20 

25 

Rotation  Rate  -  deg/sec 

60 

65 

70 

Final  Results 

Position  -  m 

-6378132.00 

1.88 

-1.25 

Velocity  -  m/sec 

-10.00 

-7.50 

-5.00 

Orientation  -  deg 

64.32 

14.84 

72.50 

Table  2.  Test  of  DR  algorithm  4,  Delta  time  is  0.5  sec 
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Initial  conditions  x  y  z 

Position -m  6378137  5  10 

Velocity  -  m/sec  30  35  40 

Acceleration  -  m/sec*sec  -80  -85  -90 

Orientation  -  deg  15  20  25 

Rotation  Rate  -  deg/sec  -60  -65  -70 

Final  Results 

Position -m  6378144.0  10.46  18.32 

3 _ 

Velocity  -  m/sec  -10.00  -7.50  -5.00 

Orientation  -  deg  -23.07  -8.68  -10.86 

Table  3.  Test  of  DR  algorithm  8,  delta  time  is  0.5  sec 

Following  Table  4  is  a  time  performance  data  comparison  between  this  thesis’s 
code  implementation  in  C-I-+  and  ARPA’s  study  code  implementation  in  C.  Those 
implementations  were  run  on  different  platforms  to  test  the  consistency  and  time  usage  of 
the  code. 


DRM 

ARPA’s  study 
SGI  Indigo  2 
lOOMHz 

ARPA’s  study 
SGI  Indigo  2 
150MHz 

Thesis  Thesis  Thesis 

Implementation  Implementation  Implementation 
SGI  Indigo  2  SGI  Indigo  2  486DLC 

150MHz  150MHz  33MHz 

Optimized 

1 

1.9  usee 

1.3  usee 

0.08  usee 

0.08  usee 

7  usee 

2 

3.8  usee 

2.6  usee 

0.96  usee 

0.40  usee 

45  usee 

3 

55  usee 

30  usee 

35.20  usee 

26.32  usee 

900  usee 

D 

57  usee 

30  usee 

36.64  usee 

26.69  usee 

980  usee 

5.9  usee 

4.0  usee 

2.08  usee 

1.12  usee 

7  usee 

6 

20  usee 

1 1  usee 

0.24  usee 

0.16  usee 

20  usee 

7 

63  usee 

33  usee 

51.60  usee 

34.56  usee 

1400  usee 

8 

71  usee 

38  usee 

56.84  usee 

44.20  usee 

1900  usee 

9 

22  usee 

13  usee 

0.88  usee 

0.40  usee 

50  usee 

Tabled.  Performance  Data 
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Observe  that  the  trend  of  the  times  spent  on  computation  were  similar  in  all  cases 
and  different  platforms.  Without  optimization  in  compiling,  this  thesis  implementation  ran 
fast  in  the  simple  cases  of  algorithms  1, 2,  5, 6,  and  9,  but  slow  in  the  complicated  cases  of 
3,  4,  7,  and  8.  Using  compilation  optimization  (with  the  flag  “-03  -mips2”),  and 
considering  the  lower  overhead  resulting  from  putting  each  class  into  its  own  file,  most  of 
the  cases  ran  faster  than  the  ARPA’s  study  code.  However,  in  the  complex  case  of 
algorithm  8,  our  code  ran  a  little  slower  due  to  redundant  matrix  multiplication  in 
computing  angular  position.  The  speed  of  the  486DLC  33Mhz  CPU  is  about  three  times 
slower  than  the  Intel  486DX2  66MHz  CPU. 

In  summation,  the  code  for  each  of  the  DR  algorithms  works  as  designed,  and  has 
been  implemented  in  the  current  NPSNET.  NPSNET  thus  takes  an  important  step  forward 
in  DIS-compliance.  NPSNET  can  move  networked  entities  from  disparate  simulators 
following  any  of  the  DR  algorithms  those  entities  specify. 
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IV.  MUNITIONS 


A.  OVERVIEW 

One  of  the  goals  of  simulation  is  to  make  the  virtual  world  resemble  the  real  world 
as  closely  as  possible.  This  chapter  describes  an  extensive  munitions  class,  built  so  that 
NPSNET  programmers  can  enhance  and  improve  the  physically-based  representation  of 
munitions  in  NPSNET.  The  class  incorporates  the  mathematics  that  will  allow  various 
munitions  to  follow  realistic  trajectories  in  the  simulated  world. 

The  kinematic  equations  v  =  v^  + at  and  x  =  x^  +  vj  +  0.5at^  calculate  velocity 

and  position  for  most  of  the  munitions,  thus  forming  the  core  of  the  base  class.  Orientation 
of  the  munition  is  derived  from  its  current  velocity  vector.  Implementation  of  the 
subclasses  of  the  base  class  is  straightforward  from  those  equations,  except  for  the  rocket 
model  which  involves  a  second  order  differential  equation. 

B.  DESIGN 

The  munition  class  is  the  base  class  containing  all  the  common  variables  and 
functions  for  its  subclasses,  see  Figure  6.  The  equations  that  govern  the  velocity  and 
position  of  the  bullet,  bomb,  and  ballistic  subclasses  are  much  the  same.  In  the  artillery 
class,  though,  we  find  the  artillery’s  necessary  initial  velocity  first,  then  it  uses  the  same 
equation  as  its  parent  class— bomb/ballistic— to  find  its  velocity  and  position  at  a  later  time. 
The  rocket  subclass  is  totally  different  from  its  peers.  Its  velocity  and  position  are  governed 
by  a  different  set  of  complex  equations.  The  approach  to  compute  the  velocity  and  position 
of  the  ballistic  missile  class  is  the  same  as  for  the  artillery  class.  That  is,  its  initial  conditions 
are  determined  before  using  the  rocket’s  velocity  and  position  functions  to  find  its  velocity 
and  position  at  a  later  time. 

The  guided  missile  seeks  and  destroys  its  target.  It  is  a  subclass  of  the  base  class. 
Incorporated  with  the  engine  class,  its  trajectory  can  be  further  refined  by  factoring  in  the 
power  curve  as  it  is  affected  by  user  input.  As  examples,  in  one  case  the  missile  may  couple 
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Figure  6:  Munition  class  structure 


with  a  constant  speed  curve  to  chase  after  the  target,  while  in  another  it  couples  with  the 
quadratic  engine  speed  equation. 

C.  IMPLEMENTATION 
1.  Orientation 

The  orientation  of  the  munition  is  converted  from  its  velocity  vector,  a  function 
residing  in  the  munition  base  class.  The  function  takes  a  velocity  vector  and  returns  angular 
position  in  degrees,  the  values  commonly  used  by  graphics  software.  This  function  checks 
all  eight  cases  of  given  (x,  y,  z)  values,  and  then  used  inverse  tangent  or  direct  assignment 
to  assign  appropriate  angular  position  for  each  direction. 

For  example,  if  the  input  vector  is  (x,  y,  z)  then: 


0  =  r  to  d  •  atan- 

Eq21 

JC 

y 

V 

0  =  r  to_d  •  atan- 

Eq22 

X 

1 

0  =  r  to  d  •  atan- 
1  X 

Eq23 

Where  r_to_d  is  the  factor  to  change  from  radian  to  degrees.  If  the  input  vector  is 


0  =  90 

X 


Eq24 


0  =  r  to  d  •  atan  -  Eq  25 

3^  “  “  z 

0  =  r  to  d  •  atan  -  Eq  26 

z  -  -  j 

2.  Bullet 

A  bullet  is  a  small  projectile  fired  from  small  arms.  To  simulate  its  path,  we  use 
equations  Eq  27  through  Eq  32.  The  gravity  constant  is  set  to  one  or  zero  to  implement  the 
bullet  with  or  without  the  effect  of  gravity,  see  Figure  7  and  equations  Eq  27  through  Eq  32. 


Figure  7:  Bullet  path 


Inputs: 

v-  initial  velocity  in  the  respective  direction. 
x-  initial  position  in  the  respective  direction. 

Outputs: 

Vj-  (0  velocity  at  time  t  in  the  respective  direction. 

JC|-  (0  position  at  time  t  in  the  respective  direction. 

cOj-  (?)  angular  position  at  time  t  in  the  respective  direction. 

The  following  equations  Eq  27  through  Eq  29  describe  the  velocity  of  the  bullet. 
Velocity  in  the  x  and  y  axes  of  the  horizontal  plane  will  not  change  from  the  initial  state 
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because  there  is  no  significant  force  acting  on  it.  We  omit  the  air  drag  and  wind  effects, 
which  are  negligible.  Gravitational  force  acts  upon  the  z  direction. 

V  =  V  Eq27 


V 

z 


V 


V 

z 

o 


yo 

-  gt  ■  gravity 


Eq28 

Eq29 


The  following  equations  Eq  30  through  Eq  32  describe  the  position  of  the  bullet. 

Position,  at  any  time,  is  based  on  initial  position  and  initial  velocity.  In  the  z  direction,  the 

bullet  drops  due  to  gravitational  force. 

X  =  X  +v  t  Eq  30 

XXX 

o  o 

X  =  X  +V  t  Eq31 

y  y  y 

^  o  ^  o 

X  -X  +v  t-O.SgP' ■  gravity  Eq32 

^  ^  o  o 


3.  Bomb 

Typically,  a  bomb  is  dropped  from  a  flying  object.  Due  to  gravity,  it  descends  from 
its  initial  position  and  accelerates  toward  the  surface  of  the  earth.  Its  path,  normally,  is  a 
parabola.  Figure  8.  The  equations  governing  its  velocity  and  position  are  much  the  same  as 
the  equations  for  the  bullet,  Eq  27  through  Eq  32.  Inputs  and  outputs  are  the  same  as  for  the 
bullet  class  above. 


Together  with  Eq  21  and  Eq  22,  the  following  equation  describes  the  velocity  of  the 

bomb. 
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V  = 
z 


Eq33 


Together  with  Eq  30and  Eq  31,  the  following  equation  describes  the  position  of  the 

bomb. 


X  =  x  +v  t-0.5gt^ 


Eq34 


These  equations  are  the  same  as  Eq  27  through  Eq  32.  The  difference  is  the 
“gravity”  constant  of  Eq  29  and  Eq  32  is  omitted,  since  bombs  are  always  dropped  from 
high  above  the  surface  where  gravitational  force  always  applies. 

4.  Ballistic 

Ballistic  behavior  applies  to  projectiles,  like  an  artillery  round,  which  are  shot  not 
directly  at  the  target,  but  along  a  parabola  whose  path  is  calculated  to  culminate  at  the  target 
location.  Figure  9.  Inputs  and  outputs  are  the  same  as  for  the  bullet  and  bomb  classes  above. 


25 


initial  position  final  position 

Figure  9:  Path  of  the  ballistic 
Equations  that  govern  the  behavior  of  ballistic  projectiles  are  much  the  same  as 
equations  for  the  bomb.  The  only  difference  between  the  two  classes  are  the  initial 
conditions,  that  is,  the  ballistic  projectile  normally  originates  from  the  ground  and  the  bomb 
from  the  air. 

a.  ArtUlery 

The  path  of  the  artillery  round  is  the  same  as  for  the  general  ballistic  round. 
The  difference  is  in  the  given  conditions.  For  the  ballistic,  we  are  given  the  initial  velocity, 
and  for  the  artillery  we  are  given  the  artillery  and  the  target  positions. 

Inputs; 

Pl(xl,  yl,  zl)  initial  position  of  the  artillery. 

P2(x2,  yl,  zl)  initial  position  of  the  target. 

Outputs:  the  same  as  for  the  bullet. 

In  the  following  case,  we  only  know  two  positions:  PI  and  P2.  We  do  not 
know  the  artillery  round’s  initial  velocity,  see  Figure  8.  If  we  did  know  the  artillery’s  initial 
velocity,  the  way  to  find  the  position  at  time  t,  P(t),  would  be  the  same  as  for  the  bomb  or 
ballistic.  So  the  approach  is  to  find  the  initial  velocity  first.  Once  we  have  the  initial 
condition,  the  equations  that  govern  the  behavior  of  the  velocity  and  the  position  at  a  later 
time  are  the  same  as  the  equations  for  the  bomb  or  ballistic. 


Figure  10:  Artillery’s  round  path 


The  following  equations  Eq  35  through  Eq  37  are  used  to  find  the  required 
nozzle  velocity  of  the  artillery  round. 

V  =  -  Eq35 

^0  t 

y-yo 

\  =  —  Hq36 

|z  — z  +0.5gt^ 

V  =  - - Eq37 

t 

In  Eq  35  through  Eq  37,  the  unknown  variable  on  the  right  is  the  time,  t,  so 
the  programmer  can  code  the  implementation  to  query  the  user  for  that  value.  In  Eq  37,  the 
top  right  part  takes  the  absolute  value  to  ensure  that  the  velocity  in  the  z  direction  is  positive 
(going  up). 

5.  Rocket 

A  rocket  is  a  type  of  munition  that  burns  fuel,  shooting  it  rearwards  to  gain  velocity 
and  altitude,  hence  distance  (Figure  11).  At  a  certain  bum-out  point,  the  rocket  will  have 
expended  its  fuel.  From  that  point,  it  enters  free-fall  as  a  ballistic  projectile.  Its  subsequent 
path  is  based  on  the  velocity  and  position  gained  in  the  bum  phase.  Therefore,  there  are  two 
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set  of  equations  governing  the  rocket’s  velocity  and  position,  the  first  applying  to  the  burn 
phase,  the  second  to  the  free-fall  phase. 

Inputs:  all  the  variables  on  the  right  side  of  the  following  velocity  and  position 
equations. 

Outputs;  the  same  as  before;  velocity,  position,  and  angular  position. 

Equation  38  states  that  the  rate  of  change  of  the  momentum  (P  =  mv)  is  equal  to  the 
sum  of  the  external  forces  that  act  on  the  object. 

I?  =  JJ’e,, 

Based  on  Eq  38,  we  can  derive  the  equation  for  the  motion  of  an  object  with  lost 
mass,  Eq  39,  where  m  is  mass,  V  is  velocity,  u  is  the  velocity  of  the  fuel,  t  is  time,  and  F  is 
the  external  force.  Eq  39  states  that  mass  multiplied  by  the  rate  of  change  of  velocity,  minus 
velocity  multiplied  by  the  rate  of  change  of  mass,  equals  the  external  force  that  applies  to 
the  object. 
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-  t 

dt  ext 


From  Eq  39,  and  if  we  assume  that  ^  =  r  where  r,  rate  of  lost  mass,  is  a  constant, 

dt 

we  can  derive  the  velocity  and  position  equations.  The  derivation  can  be  found  in  many 
math  and  physics  books  [MART84][HALL78][BOYC77]. 

Eq  40  through  Eq  42  following  are  the  velocity  equations  for  the  rocket. 


v-v  =Mln  1 - -Qt 

z  1  z  m 

O  \  o 


v-v  =  M  In  1-  — 
XXX  m 

o  I  o 


v-v  =  M  In  1-  — 
y  m 

o  \  o 


Eq  43  through  Eq  45  following  are  the  equations  for  the  position  of  the  rocket. 
u  m  (  /  .  ^  _  /  ,  N-'i 


X  -X  =  V  t- 
1  z  1 
o  o 


,  f,  rt  V,  ,  (,  rOll  i  „  „ 

o  I  V  oJi  V  oJJJ 


X  -X 

u  m 

^  X  o 
=  V  t - 

X  X 

X  r 

o 

0 

X  -X 

u  m 

=  v  t  y  ^ 

y  y 

''o 

yo  ^ 

Where: 

-  [l-ln  1-—  f 
m 

oJL  V  oJJ 


1-—  fl-ln  1-—  ’ 

m  m 

0)1  oU 


=  respective  velocity  in  each  direction 
=  respective  position  in  each  direction 


il.  =  respective  fuel  velocity  in  each  direction 
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m  =  total  mass  of  the  rocket 
o 

r  =  rate  of  lost  mass 

g  =  gravitational  constant,  9.8  lm/5 

^ .  =  respective  initial  velocity  and  position 

t  =  time 


If  we  assume  that  fuel  mass  will  burn-out  at  a  certain  time,  as  in  Eq  46,  we  can 
derive  the  bum-out  velocity  and  the  bum-out  position  as  in  equations  Eq  47  through  Eq  52. 

nij.  ,  =  m  -m  =  rt  Eq46 

fuel  o  rocket 

Where: 

nij:  ,  =  fuel  mass 
fuel 

m  ,  ^  =  mass  of  the  bare  rocket 
rocket 


The  following  equations  Eq  47  through  Eq  49  give  the  burn-out  velocity  of  the 


rocket. 


'  mA  m. 

V  -  V  =  -M  In  1  +  —  -g— 
z  z  m  °  r 

o  \  r/ 


v-v  =-Mlnl-t-  — 
X  X  1  m 

o  V  r 


V  _  V  =  -M  In  1  4-  — i 

J  y  [ 


The  following  equations  Eq  50  through  Eq  52  give  the  burn-out  position  of  the 


rocket. 
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=  V  t- 
z 

o 


1- 


1-ln 


Eq50 


JC  -X 

z  z 


u  m 
z  o 


( 

\^oJ 


(  mM 
1 


V 


m 


r)l 


2^ 


( 

J 


u  m 

*  ^  i  1 

X  -X  =  V  t - 1 1  - 

XX  X  r  \ 

o  o 


v%y 


1-ln 


m. 


Y* 

1+^ 


m 


rJ} 


Eq51 


X  -X 

y 


V  t- 


u  m 
yo 


1- 


( m  \ 

r 

.m 
V  oy\_ 


1-ln 


m. 


f 

1  +  ^ 


m 


rJ\ 


Eq52 


a.  Ballistic  missile 

The  inter-continental  ballistic  missile  (ICBM)  is  a  type  of  weapon  which 
uses  rocket  thrust  and  has  the  capability  of  flying  directly  to  a  given  target  position.  Figure 
12.  The  equations  governing  its  velocity  and  position  are  the  same  as  for  the  artillery  class. 
We  determine  the  ICBM  initial  conditions.  Then  we  use  the  parent  class  velocity  and 
position  equations  to  find  the  velocity  and  position  at  a  later  time. 


Figure  12:  Rocket  shot  at  predetermined  target. 
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6.  Guided  missile 


A  guided  missile  is  a  projectile  with  the  capability  of  chasing  and  homing 
in  on  a  moving  target,  Figure  13.  Our  general  solution  for  the  moving  target  will  apply  to 
the  fixed  target  case,  too.  In  the  moving  target  case,  the  missile  must  change  its  velocity 
vector  every  delta  time.  At  =  t2  -  tl.  At  the  end  of  each  delta  time,  we  know  the  position 
and  velocity  of  both  the  missile  and  target.  We  use  that  as  the  initial  conditions  for  the  next 
time  interval.  We  repeat  the  process  until  the  missile  is  within  a  detonation  range  from  the 
target. 

The  following  equation  Eq  53  computes  the  new  missile  velocity  for  every  At. 


missile 


P  —P 

target  missile  „  , 

p  _p - 

target  missilel 


Eq53 


Where: 


'^tmget  =  l»g«P“Won 


^  .  .,  =  new  missile  velocity 

missile 

^missile  =  missile  position 

Speed  =  speed  of  the  missile  from  previous  delta  time 
Following  is  the  position  equation  for  the  guided  missile 

X.  =  X.  +v.  At  +  0.5a.{At)^  Eq54 

I  Iq  ' 

Where  i  =  is  the  respective  direction,  x  ■  =  next  position,  x.  =  current  position 

I  iQ 

direction,  v .  =  velocity  at  current  position,  a  ■  =  acceleration,  and  At  =  given  delta  time. 

Iq  » 

Remember  that  v .  is  recalculated  every  At . 

*0 

D.  RESULTS 

As  in  Chapter  3,  each  case  of  the  munitions  class  is  coded  and  tested  separately.  The 
code  for  each  case  resides  in  separate  header  and  body  files.  This  divide-and-conquer 
strategy  facilitated  the  coding  and  testing  process  tremendously.  Potential  problems  and 
errors  were  contained  within  sections  of  code  of  manageable  size.  Once  the  problems  of  the 
parent  class  are  ironed  out,  those  parent  class  attributes  inherited  by  the  child  classes  need 
not  be  validated  and  debugged  again  for  each  child  class.  With  the  C++  hierarchy  structure, 
problems  are  isolated  and  localized,  and  not  propagated. 

The  following  two  pseudo-code  procediures  outline  the  testing  of  the  munitions 
subclasses,  and  the  guided  missile  class  which  is  handled  separately: 

Following  is  the  main  function  used  to  test  the  munition  subclasses: 

- begin  main  procedure - 

1.  Instantiate  an  object 

2.  While  object  does  not  hit  the  target,  or  the  ground 

a.  Move  the  object 

b.  Get  the  object  position  (for  rendering) 

c.  Get  the  object  velocity 
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d.  Get  the  object  angular  position  (for  rendering) 

e.  Check  whether  the  object  has  hit  the  target  or  ground 

(if  it  hits,  exit  the  loop;  if  not,  continue) 

f.  Check  whether  the  object  has  missed  the  target 

(if  so,  exit  the  loop) 

■  end  while  loop - 

— end  main  procedure - 


The  following  loop  describes  the  coding  and  testing  of  the  guided  missile  class: 

- begin  main  procedure - 

1.  Set  up  initial  conditions  (delta  time,  missile  position,  target  position) 

2.  While  the  target  has  not  been  hit 

a.  Compute  missile  velocity,  angular  position 

b.  Move  the  missile  for  a  given  delta  time 

c.  Get  the  next  target  position  (from  DR) 

d.  Get  the  next  speed  (constant  or  from  speed  curve) 

e.  Determine  if  the  missile  hit  the  target  within  a  given  range 

(if  so,  exit  the  loop) 

3.  Return  to  Step  2 

- end  main  procedure - 

Note:  If  the  delta  time  is  not  small  enough  and  the  detonated  range  not  large  enough, 
the  definition  of  a  hit  required  by  step  e  will  never  be  satisfied. 


As  the  two  pseudo-code  procedures  above  describe,  the  munition’s  state  is  started 
with  initial  conditions  and  then  recalculated  after  a  predetermined  small  delta  time,  with 
checks  of  whether  the  target  has  been  hit  or  the  munition  has  impacted  the  ground.  If  one 
of  those  conditions  is  trae,  the  loop  is  exited  and  the  procedure  terminated. 

In  summation,  the  code  works  as  designed.  Each  class  returns  its  output  in  a  well- 
behaved  manner.  The  code  is  reliable,  and  the  performance  excellent,  even  on  the  low-end 
Intel  80486-based  PC  used  for  testing. 
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V.  ENGINES 


A.  OVERVIEW 

This  chapter  describes  mathematical  models  for  engine  power.  The  actual 
implementation  of  engine  models  such  as  the  gas  turbine,  motor,  and  rocket  are  not 
essential  to  the  virtual  world  of  NPSNET.  However,  it  is  important  that  the  simulation 
looks  and  feels  right  to  the  user.  Therefore,  our  goal  was  to  design  our  equations  only  to  a 
level  of  detail  that  will  have  a  beneficial  effect  on  the  simulation. 

Figure  14  shows  all  of  the  mathematical  models  which  will  be  discussed  in  this 
chapter.  Depending  on  the  model  we  use,  engine  speed  will  follow  the  mathematical  curve 
precisely,  either  the  logarithmic,  sinusoidal,  linear,  exponential,  rocket,  or  quadratic  curve. 

Note  that  we  are  only  concerned  with  increasing  or  decreasing  speed.  The  direction 
of  any  object  employing  our  engine  model  depends  on  user  input  during  the  simulation. . 


Figure  14:  Models  for  Engines 


In  Figure  14,  s  represents  the  speed  of  the  entity  and  r  is  the  range  of  the  input 
throttle.  Variable  r  represents,  for  instance,  pressing  the  gas  pedal  in  a  car.  As  r  increases, 
so  does  the  speed  of  the  car  and  vice  versa,  is  the  maximum  speed  that  the  object  can 
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attain,  and  is  the  maximum  range  number  that  the  user  will  supply.  For  example,  r 
might  run  from  0  to  100. 

B.  DESIGN 

The  engine  class  is  the  base  class.  It  contains  the  common  variables  the  subclasses, 
which  are  derived  from  the  base  class,  will  use.  Most  subclasses  have  two  cases 
implemented,  one  each  for  increasing  and  decreasing  speed.  I  and  D  in  Figure  15  shows 
that  all  subclasses  implement  both  cases,  except  for  the  rocket,  which  will  not  follow  the 
same  curve  for  decreasing  speed. 


Figure  15:  Engine  Classes 


C.  IMPLEMENTATION 

Following  are  the  symbols  we  use  to  describe  the  mathematical  model  for  engine 

speed: 

s  =  speed,  which  extends  from  0  to 
r  =  range,  which  extends  from  0  to 
a,  b,  c,  d  are  constants 
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Our  goal  is  to  find  the  simplest  equations  that  exhibit  the  characteristics  of 
increasing  and  decreasing  speed  for  each  subclass.  Those  equations  must  satisfy  the  two 
initial  conditions,  (0, 0)  and 

We  use  r  as  an  input  number  which,  when  processed  through  the  appropriate 
equation,  maps  to  a  value  in  the  speed  range  (Figure  16).  For  example,  with  equation 
s  =  2r ,  if  r  equals  2,  then  s  is  4.  In  the  code,  the  toolkit  user  chooses  engine  type  and 
specifies  and  .  After  that,  the  appropriate  s  is  returned  for  any  given  r,  and  vice 
versa. 


Following  is  the  derivation  of  possible  equations  to  represent  increasing  and 
decreasing  speed: 

1.  Increasing  and  decreasing  speed 
a.  Linear 

In  this  case,  engine  speed  follows  a  simple  linear  equation.  Consider  the 
linear  equation  for  speed 

s-a  =  slope  (r-b)  Eq55 

Where  slope  is  the  slope  of  the  linear  cmve.  In  this  case,  with  two  initial 

conditions,  (0, 0)  and  (s  ,r  )  ,  our  speed  equation  is  (See  Figure  17  and  Eq  49): 
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max 

s  =  - r 

r 

max 


S  max,  r  max 


0, 0  range 


Figure  17:  Linear  approximation 


The  range  equation  is: 


max 

r  =  - S 

^max 


b.  Quadratic 

In  this  case,  the  speed  equation  can  follow  two  quadratic  curves-the  solid 
and  dotted  lines  in  Figure  18-to  increase  or  decrease  speed. 


Figure  18:  Quadratic  approximation 


We  can  derive  the  speed  equation  for  the  solid  line  as  follows: 

With  the  general  quadratic  equation  s  =  ar^ ,  then  according  to  Figure  10 
with  our  two  initial  conditions  (0,0)  and 


Eq58 


Eq59 


We  can  also  derive  the  speed  equation  for  the  dotted  line  as  follows.The 
general  quadratic  equation  is  S  —  S  —  —a{r  —  f  ^  •  With  two  initial  conditions, 

^fYldX 

(0,0)  and  (i„„.r„„;t).thenfl  =  - - ^ .  The  speed  equation  is: 


Eq60 


Eq61 
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speed  can  increase  in  a  higher  degree  polynomial  but,  though 
straightforward,  is  not  implemented  here  owing  to  lack  of  necessity.  The  general  format  for 
coding  and  using  the  code  would  be  the  same. 

c.  Logarithmic 

In  this  case,  speed  increases  along  a  logarithmic  curve.  The  general  equation 

for  this  case  is: 

s  =  alog  (l+r)  Eq62 


Figure  19:  Logarithmic  approximation 
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r  =  e 


^log(l+r^J 

^  max 


-1 


Eq64 


d.  Exponential 

In  this  case,  speed  increases  along  an  exponential  curve.  The  general 
equation  is  s  +  b  =  with  two  initial  conditions  of  (0,  0)  and  •  We 


derive  b  =  \  and  s  +  1 
max 


nr 

e  then  a 


log  (5  +1) 

c?  V  mnx  ' 

- — - .  So  our  equation 

r 

max 


for  the  speed  is: 

s+l  =  e 

The  range  equation  is: 


Eq65 


r 


r 

max 


log  {s 


max 


+  1) 


log  (5+1) 


Eq66 


Figure  20:  Exponential  approximation 
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0  =  - 
r 


curve. 


e.  Sinusoidal 

In  this  case,  speed  follows  the  sinusoidal  curve.  See  Figure  13. 


Figure  21:  Sinusoidal  approximation 


Where  0  € 


,  we  can  map  input  to  0  easily  with  equation 


.  The  range  equation  is: 


r  = 


2 

-r 

jj-  max 


asin 


s 

p 

max 


Eq68 


/.  Rocket 

In  this  case,  the  rocket  model  in  Chapter  4  governs  the  behavior  of  the  speed 
In  general,  the  speed  equation  for  the  rocket  is: 


S  =  —Min 


^  rt^ 

— 


m 


V 


Eq69 


o) 
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Where  s  =  speed,  u  =  fuel  speed,  t  =  rate  of  burning  fuel,  r  =  range  of  input, 
=  total  mass  of  the  object,  and  rt  =  mass  of  the  fuel. 

Figure  22  describes  the  behavior  of  the  speed  curve.  The  speed  increases 
very  rapidly  to  a  very  large  value.  Initial  conditions  are  (0,0)  and  . 


We  use  Eq  62  to  implement  the  code  as  follows.  The  programmer  supplies 
as  input  u,  m  ,  r,  and  r  .  From  these  initial  values  and  by  controlling  the  u  value,  we 

assign  to  a  value  we  want.  With  known  values  of  and  the  two  equations, 

Eq  62  and  Eq  63,  we  can  now  simulate  the  rocket  engine. 

The  range  equation  is: 

m 

r  =  _  ( 1  _  Eq70 
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2.  Increasing  and  decreasing  speed  on  different  curves 

Usually,  for  the  sake  of  simplicity,  the  decreasing  speed  curve  for  an  engine  is  the 
same  as  the  increasing  case.  But  we  can  always  use  different  speed  curves  on  the  same 
engine.  One  of  the  problems  we  encountered  is  in  how  to  represent  the  input  range  of  the 
throttle  to  the  user.  If  we  use  the  two  curves  in  Figure  23,  then  for  one  value  of  si  we  have 
two  range  values,  r  1  and  r2.  Which  range  value,  then,  do  we  return  to  the  user?  One  solution 
is  to  only  return  the  r  value  of  the  increasing  speed  curve.  Then  when  engine  speed 
decreases,  we  map  the  decreasing  speed  to  the  equivalent  increasing  speed,  and  from  that 
speed  value  and  the  range  equation  of  the  increasing  case,  we  can  find  value  r. 


D.  RESULTS 

The  code  is  organized  into  a  base  class  and  dependent  subclasses.  The  base  class 
holds  the  variables  the  other  subclasses  will  use,  such  as  maximum  speed,  maximum  range, 
speed,  and  range.  Each  subclass  has  two  central  functions.  One,  get_speed(),  returns  the 
speed  for  a  given  range  value,  and  the  other,  get_range(),  returns  the  range  for  a  given 
speed. 

To  incorporate  the  code  to  use  for  each  case,  the  user  includes  a  specific  set  of  files. 
For  example,  to  use  the  linear  engine  curve,  include: 
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-  engine.cpp 

-  linear.cpp  (the  name  of  the  curve  to  use) 

-  a  main  function 

The  code  was  tested  separately  for  each  curve.  Two  loops  were  created  for  return 
values  of  speed  and  range.  For  return  speed  values,  the  range  value  extended  from  0  to  100 
in  steps  of  10.  For  range  return  values,  the  speed  value  extended  from  0  to  1000  in  steps  of 
100.  See  the  code  in  Appendix  C  for  details. 

For  the  rocket  engine  model,  the  programmer  should  be  careful  with  the  given  value 
of  fuel  speed  u,  rate  of  burning  fuel  t,  and  the  total  mass  .  Remember  that  tr^^  is  the 

total  fuel,  it  should  be  about  half  of  the  total  mass  .  The  best  value  for  maximum  range 
^ max  ^bout  100,  and  use  a  fuel  speed  value  u  to  control  the  maximum  speed  value  . 

The  return  values  for  each  case  behave  according  to  the  name  of  the  case.  For 
example,  the  linear  model  returns  values  in  both  speed  and  range  linearly,  and  the 
exponential  model  returns  speed  values  exponentially  and  its  range  values  logarithmically 
as  its  equations  describe. 

In  general,  all  models  worked  very  well.  The  particular  case  is  the  rocket  engine. 
With  its  physically-based  model,  it  should  be  used  for  all  jet  engines  in  the  simulation  to 
enhance  the  realistic  effects  of  our  simulated  world. 

With  the  division  of  separate  subclasses  branching  from  the  base  class,  the  usage, 
testing,  and  enhancement  of  the  code  is  straightforward.  In  sununation,  the  models  work 
very  well  in  our  simulated  world. 
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VI.  RESULTS  AND  CONCLUSIONS 


A.  RESULTS 

The  divide-and-conquer  method  was  used  to  code  and  test  the  implementation 
described  in  this  thesis,  that  is,  each  class  was  written  in  a  separate  file,  segregated  into  a 
self-contained  chunk  of  code.  Therefore,  potential  problems  and  bugs  were  localized  and 
the  development  process  was  reduced  to  small,  manageable,  straightforward  steps.  Once  a 
particular  class  was  validated,  it  required  no  further  examination.  Its  derived  classes  did  not 
retirni  further  problems  in  exchange  for  the  benefits  they  received.  Instead,  they  gained  the 
benefit  of  being  built  onto  a  solid  foundation.  In  general,  the  code  is  reliable,  performs  well, 
and  is  easy  to  use.  All  of  the  interface  functions  have  a  standardized  format,  and  follow  a 
standard  naming  convention. 

1.  Dead  Reckoning 

The  design  and  implementation  of  the  nine  DR  algorithms  in  this  thesis  are 
complete  in  the  sense  that  they  cover  all  possible  cases  in  both  world  and  body  coordinates. 
In  each  perspective  coordinate  system,  the  implemented  algorithms  cover  cases  ranging 
from  the  simplest,  a  fixed  entity,  to  the  most  detailed,  an  entity  with  rotation,  velocity,  and 
acceleration. 

Most  equations  used  to  compute  the  results  are  simple,  except  for  body  coordinate 
algorithms  7  and  8,  which  are  somewhat  long  and  complex.  However,  all  algorithms  were 
coded  with  an  efficiency  as  near  to  optimal  as  possible.  For  more  examination  of  the 
applicable  uses  of  these  algorithms,  consult  the  ARPA  study  for  more  details.  [IST93] 

2.  Munitions 

The  munitions  class  in  this  thesis  formed  a  basis  for  most  of  the  munitions.  It 
accounts  for  the  bullet,  bomb,  ballistic,  artillery,  rocket,  and  guided  missile.  Further,  the 
user  can  combine  these  munitions  to  form  different  types  of  weapon  with  each  distinct 
munition  performing  a  single  phase  of  some  composite  weapon.  For  example,  the  bomb 
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model  could  be  used  to  drop  a  missile  from  an  airplane,  and  then  when  it  reaches  a  certain 
altitude,  the  guided  missile  model  could  take  over  to  govern  a  final  stage  wherein  the 
missile  seeks  a  target. 

The  guided  missile  mode  of  operation  can  be  tail-chase  or  predicted-intercept 
depending  on  what  kind  of  target  position  the  user  provides  to  compute  the  missile  velocity 
vector.  K  the  user  provides  the  current  position  of  the  target,  then  the  missile  mode  will  be 
tail-chase,  because  the  target  will  have  moved  to  the  next  position  at  the  next  delta  time.  If 
the  user  provides  the  next  dead  reckoned  position  of  the  target,  the  missile  mode  will  be  the 
predicted-intercept  missile. 

Furthermore,  the  guided  missile  class  can  couple  with  any  of  the  engine  models  in 
Chapter  5.  Recall  that  the  direction  of  the  missile  is  based  on  the  position  of  the  missile  and 
the  position  of  the  target.  If  the  speed  from  the  engine  equation  for  every  delta  time  is 
supplied  to  the  munition,  then  the  missile  will  chase  after  the  target  with  the  increasing 
speed  of  the  engine  curve. 

The  rocket  model  is  based  on  Newton’s  equation.  Its  velocity  and  position  equation 
are  final.  Together  with  the  bomb  equations,  the  rocket  equations  will  form  the  basis  for  the 
ballistic  missile  model.  However,  the  ballistic  missile  model  was  not  implemented  due  to 
time  constraints  and  its  complexity. 

3.  Engine 

Except  for  the  rocket  model  which  is  physically-based,  the  remaining  engine 
models  are  mathematically-based.  But  thorough  testing  of  all  the  models  reveals  negligible 
difference  between  physically-  and  mathematically-based  models. 

In  physically-based  models,  we  use  Newton’s  force  law  to  get  a  mathematical 
expression  of  the  speed  curve,  whereas  in  mathematical  models,  speed  curves  are  encoded 
based  on  popular  mathematical  equations.  In  both  cases,  the  initial  condition  (0,  0),  the 
maximum  speed,  and  the  maximum  range  are  used  as  boundary  values.  Our  model  becomes 
a  boundary  value  problem,  that  is,  the  fitting  of  a  curve  between  two  points. 
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In  testing,  since  the  quadratic  or  exponential  curve  behaves  very  much  the  same  as 
the  logarithnruc  rocket  curve,  we  can  use  them  to  substitute  for  the  rocket  speed  curve. 
Overall,  our  models  behaved  realistically  when  tested. 

B.  CONCLUSIONS 

The  code  described  in  this  thesis  keeps  each  separate  class  independent  from  each 
other.  It  is  short  and  the  performance  fast.  The  names  of  the  functions  are  standardized, 
making  the  toolkit  easy  to  use.  The  toolkit  was  tested  thoroughly  to  ensure  reliability. 

The  work  of  this  thesis  cannot  be  considered  groundbreaking,  though.  But  it  does 
craft  existing  knowledge  into  an  organized  and  coherent  piece  of  work.  Its  design  and 
implementation  are  well  thought  out,  and  its  performance  in  the  current  system  is  excellent. 
Furthermore,  it  lays  a  solid  foundation  for  further  study. 

C.  SUGGESTIONS  FOR  FUTURE  WORK 

After  building  the  toolkit,  the  advantage  of  concise,  independent  functions  that  can 
be  coupled  with  any  entity  models,  was  clear.  We  suggest  that  the  work  for  suspension, 
control,  steering  models,  and  other  dynamic  entity  attributes,  be  implemented  in  a  manner 
consistent  with  the  toolkit. 
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